home *** CD-ROM | disk | FTP | other *** search
/ Shareware Grab Bag / Shareware Grab Bag.iso / 090 / byt1186b.arc / APPROX.LBR / CHBEVL.C next >
Text File  |  1986-04-11  |  2KB  |  83 lines
  1. /*                            chbevl.c
  2.  *
  3.  *    Evaluate Chebyshev series
  4.  *
  5.  *
  6.  *
  7.  * SYNOPSIS:
  8.  *
  9.  * int N;
  10.  * double x, y, coef[N], chebevl();
  11.  *
  12.  * y = chbevl( x, coef, N );
  13.  *
  14.  *
  15.  *
  16.  * DESCRIPTION:
  17.  *
  18.  * Evaluates the series
  19.  *
  20.  *        N-1
  21.  *         - '
  22.  *  y  =   >   coef[i] T (x/2)
  23.  *         -            i
  24.  *        i=0
  25.  *
  26.  * of Chebyshev polynomials Ti at argument x/2.
  27.  *
  28.  * Coefficients are stored in reverse order, i.e. the zero
  29.  * order term is last in the array.  Note N is the number of
  30.  * coefficients, not the order.
  31.  *
  32.  * If coefficients are for the interval a to b, x must
  33.  * have been transformed to x -> 2(2x - b - a)/(b-a) before
  34.  * entering the routine.  This maps x from (a, b) to (-1, 1),
  35.  * over which the Chebyshev polynomials are defined.
  36.  *
  37.  * If the coefficients are for the inverted interval, in
  38.  * which (a, b) is mapped to (1/b, 1/a), the transformation
  39.  * required is x -> 2(2ab/x - b - a)/(b-a).  If b is infinity,
  40.  * this becomes x -> 4a/x - 1.
  41.  *
  42.  *
  43.  *
  44.  * SPEED:
  45.  *
  46.  * Taking advantage of the recurrence properties of the
  47.  * Chebyshev polynomials, the routine requires one more
  48.  * addition per loop than evaluating a nested polynomial of
  49.  * the same degree.
  50.  *
  51.  */
  52. /*                            chbevl.c    */
  53.  
  54.  
  55. /* Cephes Math Library Release 1.1:  March, 1985
  56.  * Copyright 1985 by Stephen L. Moshier
  57.  * Contributed to BIX for personal, noncommercial use only.
  58.  * Direct inquiries to 30 Frost Street, Cambridge, MA 02140 */
  59.  
  60. double chbevl( x, array, n )
  61. double x;
  62. double array[];
  63. int n;
  64. {
  65. double b0, b1, b2, *p;
  66. int i;
  67.  
  68. p = array;
  69. b0 = *p++;
  70. b1 = 0.0;
  71. i = n - 1;
  72.  
  73. do
  74.     {
  75.     b2 = b1;
  76.     b1 = b0;
  77.     b0 = x * b1  -  b2  + *p++;
  78.     }
  79. while( --i );
  80.  
  81. return( 0.5*(b0-b2) );
  82. }
  83.